#include "imageprocessor.h"

ImageProcessor::ImageProcessor()
{
    height = width = 0;

    // Default output file
    outFile = QString("./output/seed_germination.csv");
    // Default dimension value for GLCM vector
    dim = DISTLEVEL*FEAT;
}

ImageProcessor::~ImageProcessor()
{

}

void ImageProcessor::deleteMem()
{
    for(int i = 0; i < nPts; i++)
        delete [] dataPts[i];
    delete [] dataPts;
}

/*
 *  In 'v' image, apply wavelet transform for each region containing one seed. Then calcaulate
 *  coccurance texture features on suband 4, 5, 6 based on integer DWT coefficients.
 */
void ImageProcessor::analyzeWaveletRegion(const Mat& v, const QList<SeedNode *>& regions, bool runNNFlag)
{
    if (regions.isEmpty()) return;
    int count = regions.size();

    SeedNode *node = 0;

    // Create one dimensional array of wsize x wsize that
    // stores wavelet coefficients value for each seed region
    float *f_data = 0;
    // Integer type of wavelet coefficients
    int *i_data = 0;

    // Subbands of wavelet coefficients, and calculate FEAT number of features
    // for each subband, then repeat these steps for all distances [0..DISTLEVEL)
    // Therefore total number of dimensions is NUMBAND*FEAT.
    allocPts(count);

    width = v.cols;
    height = v.rows;

    for(int i = 0; i < count; i++)
    {
        // Each top-left corner point of seed region
        node = (SeedNode *)regions.at(i);
        // Make sure node points to a valid seednode
        if(node == 0) continue;

        f_data = new float[A_SIZE];
        i_data = new int[A_SIZE];

        // extract data from each region within image and store them in
        // one dimensional float array 'f_data', then perform wavelet
        // transform on coefficient, convert obtained float coefficients
        // to integer type. Then calculate coccurance GLCM
        getRegionData(f_data, v.data, v.step, node->pos);
        waveletTransform(f_data, i_data);
        waveletCoccurance(i_data, dataPts[i]);

        delete [] i_data;
        delete [] f_data;
        i_data = 0;
        f_data = 0;
    }

    // output dataPts for training
    if(!runNNFlag) outputSummary(regions);
}

/*
 * 'v' is of gray format, usually in VALUE channel of HSV color space. All pixel value will be quantized to
 * 'LEVEL' and 'v' has regions enclosing seeds ready to be analysed, identified by 'regions'
 * In this method, an image of seeds marked in regions need to process. All regions are squares
 * of size WSIZE x WSIZE. Each region is of type seednode that contains information related to
 * to each seed such as location and germination result: 0 for non-germination and 1 for germination
 *
 */
void ImageProcessor::analyzeRegion(const Mat& v, const QList<SeedNode *>& regions, bool runNNFlag)
{
    if (regions.isEmpty()) return;
    int count = regions.size();

    SeedNode *node = 0;

    // Create one dimensional integer array of size x size, which
    // stores pixel value for each seed region
    int *data = 0;
    // Create a coccurance for each array.
    // Coccurance(h, w, n, level). In this case, size = 64
    // Coccurance (height, width, ng);
    Coccurance *cc = new Coccurance(WSIZE, WSIZE, LEVEL);

    // 'FEAT' stands for number of dimension of GLCM
    // 'count' stands for number of data points. Repeat
    // the allocation for all distances [0..DISTLEVEL)
    // Therefore, total number of features will be FEAT*DISTLEVEL
    allocPts(count);

    width = v.cols;
    height = v.rows;

    // Iterate each region, extract gray level pixel values,and calculate texture feature vector
    // then save calculated results back to data points
    for(int i = 0; i < count; i++)
    {
        // Each top-left corner point of seed region
        node = (SeedNode *)regions.at(i);
        // Make sure node points to a valid seednode
        if(node == 0) continue;

        data = new int[A_SIZE];

        // extract data from iplimage and store them in integer array 'data'
        getRegionData(data, v.data, v.step, node->pos);
        distCoccurance(data, dataPts[i], cc, 0);

        delete [] data;
        data = 0;
    }

    delete cc;
    // output dataPts for training
    if(!runNNFlag) outputSummary(regions);
}

/*
 * Get data from 'src' and store into float array 'dst'
 */
void ImageProcessor::getRegionData(float *dst, unsigned char *src, int widthstep, CvPoint &p)
{
    int i, j;

    int x = p.x;
    int y = p.y;

    int y_end = y + WSIZE;
    int x_end = x + WSIZE;

    float value = 0.0;

    for(i = y; i < y_end; i++)
    {
        for(j = x; j < x_end; j++)
        {
            // If region is out of boundary, assign 0's
            // Otherwise assign image value to it.
            if(i >= height || j >= width || i < 0 || j < 0)
                value = 0.0;
            else
                value = (float) src[i*widthstep + j];

            dst[(i-y)*WSIZE+j-x] = value;
        }
    }
}
/*
 * Get data from 'src' and store into int array 'dst'
 */
void ImageProcessor::getRegionData(int *dst, unsigned char *src, int widthstep, CvPoint &p)
{
    int i, j;
    int value = 0;
    int x = p.x;
    int y = p.y;

    int y_end = y + WSIZE;
    int x_end = x + WSIZE;

    for(i = y; i < y_end; i++)
    {
        for(j = x; j < x_end; j++)
        {
            // If region is out of boundary, assign 0's
            // Otherwise assign image value to it.
            if(i >= height || j >= width || i < 0 || j < 0)
                value = 0;
            else
                value = (int) src[i*widthstep + j];

            dst[(i-y)*WSIZE+j-x] = value;
        }
    }
}

/*
 * idata is integer type coefficients of wavelet transform
 * of size (w * h).
 */
void ImageProcessor::waveletCoccurance(int *idata, double* p)
{
    // Divider of wavelet subbands
    int div = 2 * WLEVEL;

    // Depending on levels of wavelet transform, coccurance height
    // and width will be same as size of each wavelet subband
    // Coccurance(int h, int w, int level);
    Coccurance *cc = new Coccurance(WSIZE/div, WSIZE/div, LEVEL);

    // 4th suband
    ////////////////////////////////////////////
    distCoccurance(idata+WSIZE/div, p, cc, 0);

    // 5th suband
    // distCoccurance(idata+A_WSIZE/div+WSIZE/div, p, cc, 1);

    // 6th suband
    ////////////////////////////////////////////
    distCoccurance(idata+A_SIZE/div, p, cc, 1);

    delete cc;
}

// Calculate GLCM coccurance for integer data and store obtained GLCM
// features to datapts 'pts'. Based on number of distance, 'DISTLEVEL'
// repeat above steps for each distance value. 'group' determine which
// group of distance GLCM features will be calculated. For non-wavelet
// features, group [0..3] stands for GLCM features obtained at each
// distance value. For wavelet features, group [0..3] for band 4 and group
// [4..7] for band 6
void ImageProcessor::distCoccurance(int *data, double *pts, Coccurance *cc, int g)
{
    // Distance value is based on pow(2, i) where i
    // loops within [0..DISTLEVEL)
    int d = 0;
    int gp = g * DISTLEVEL;

    // Loop through each distance value
    for(int i = 0; i < DISTLEVEL; i++)
    {
        d = (int)pow(2., i);

        cc->setParam(data, d);
        fillPts_Light(cc, pts, gp+i);

        cc->reset();
    }
}

// Do wavelet transform for each region of image of
// size (WSIZE x WSIZE)
void ImageProcessor::waveletTransform(float *data, int *idata)
{
    // Wavelet(int h = 0, int w = 0, int l = 1);
    Wavelet *wave = new Wavelet(WSIZE, WSIZE, WLEVEL);
    // Perform forward wavelet transform on data
    wave->transform(data);
    // Convert float coefficients to int type
    convertToInt(data, idata);
    // delete 'wave' and clear memory
    delete wave;
}

// Allocate dataPts with dimension 'd' and number of
// data points 'n'. For gray level coccurance, number
// of dimension is distance level multiply with number
// GLCM features. For wavelet coccurance, number of
// dimension is distance level multipl with number of
// GLCM features as well as number of suband chosen
void ImageProcessor::allocPts(int n)
{
    nPts = n;

    dataPts = new double* [n];
    for(int i = 0; i < n; i++)
        dataPts[i] = new double[dim];
}

/********************************************************************************
 ********************************************************************************
 *
 *                    Utilities member functions
 *
 ********************************************************************************
 ********************************************************************************
 */
// --------------------------------------------------
// Fill up data points from calculated coccurance
// --------------------------------------------------
void ImageProcessor::fillPts_Light(Coccurance* cocc, double* p, int index)
{
    int start = index * FEAT;
    int end = (index+1) * FEAT;

    for(int i = start; i < end; i++)
    {
        // Get the average statistics of all 4 angles/directions
        p[i] = cocc->avg_stat(i%FEAT);
    }
}

/*
 *  Output data of training set for neural network.
 *
 */
void ImageProcessor::outputSummary(const QList<SeedNode *>& loc)
{
    try
    {
        ofstream fout;
        fout.open(outFile.toStdString().c_str(), ofstream::app);

        for(int i = 0; i < nPts; i++)
        {
            double* p = dataPts[i];

            for(int j = 0; j < dim; j++)
            {
                fout << p[j] << ",";
            }
            fout << loc.at(i)->val;
            fout << "\n";
        }

        fout.close();
    }
    catch(exception& e)
    {
        qDebug() << e.what();
    }
}

// Convert 'WLEVEL' wavelet coeffecients to integer of range [0, LEVEL-1]
// Both 'coe_int' and 'coe' are same dimension
void ImageProcessor::convertToInt(float *coe, int *coe_int)
{
    float max4=0.0, min4=0.0, max5=0.0, min5=0.0, max6=0.0, min6=0.0;    // min-max value
    float dv4=0.0, dv5=0.0, dv6=0.0;  // delta
    float t = 0.0;
    int p = 0;

    // divider of image size depending on level of wavelet transform
    int div = 2 * WLEVEL;

    int divEnd = WSIZE / div;
    int levEnd = WSIZE / WLEVEL;

    // findMinMax(float *data, float &max, float &min, int start_h, int start_w, int endw, int endh)
    findMinMax(coe, max4, min4, 0, divEnd, levEnd, divEnd);
    findMinMax(coe, max5, min5, divEnd, divEnd, levEnd, levEnd);
    findMinMax(coe, max6, min6, divEnd, 0, divEnd, levEnd);

    dv4 = (max4 - min4)/(float)LEVEL;
    dv5 = (max5 - min5)/(float)LEVEL;
    dv6 = (max6 - min6)/(float)LEVEL;

    // Upper left subband
    for(int i = 0; i < divEnd; i++)
    {
        for(int j = 0; j < divEnd; j++)
        {
            coe_int[i*WSIZE+j] = 0;
        }
    }
    //4th subband
    for(int i = 0; i < divEnd; i++)
    {
        for(int j = divEnd; j < levEnd; j++)
        {
            t = coe[i*WSIZE+j];
            // range of p [0, LEVEL-1]
            p = (int) floor(((t - min4)/dv4));
            coe_int[i*WSIZE+j] = (p==LEVEL)?(LEVEL-1):p;
        }
    }
    //5th subband
    for(int i = divEnd; i < levEnd; i++)
    {
        for(int j = divEnd; j < levEnd; j++)
        {
            t = coe[i*WSIZE+j];
            // range of p [0, LEVEL-1]
            p = (int)floor(((t - min5)/dv5));
            coe_int[i*WSIZE+j] = (p==LEVEL)?(LEVEL-1):p;
        }
    }
    //6th suband
    for(int i = divEnd; i < levEnd; i++)
    {
        for(int j = 0; j < divEnd; j++)
        {
            t = coe[i*WSIZE+j];
            // range of p [0, LEVEL-1]
            p = (int) floor(((t - min6)/dv6));
            coe_int[i*WSIZE+j] = (p==LEVEL)?(LEVEL-1):p;
        }
    }
}

void ImageProcessor::findMinMax(float *data, float &max, float &min, int start_h, int start_w, int endw, int endh)
{
    float t = 0.0;

    for(int i = start_h; i < endh; i++)
    {
        for(int j = start_w; j < endw; j++)
        {
            // Data is one dimensional array of size WSIZE x WSIZE
            t = data[i*WSIZE+j];

            if(max < t)
                max = t;
            if(min > t)
                min = t;
        }
    }
}

// 2 for wavelet and 0 for non-wavelet
void ImageProcessor::setWavelet(bool checked)
{
    if(checked)
    {
        dim = NUMBAND*DISTLEVEL*FEAT;
        outFile = "./output/wavelet_seed_germination.csv";
    }
    else
        outFile = "./output/seed_germination.csv";
}
